% Some helpful macros:
%
%\insertfigure{example_figure.png}{Figure caption goes here.}
%\ref{example_figure.png}
%
%\insertsubfigurepair
%    {label}{Shared caption}
%    {filename1.png}{subcaption1}
%    {filename2.png}{subcaption2}
%\ref{label} \ref{filename1.png} \ref{filename2.png}
%
%\insertfigurepair
%    {filename1.png}{caption1}
%    {filename2.png}{caption2}
%\ref{filename1.png} \ref{filename2.png}
%
%\inserttable{tableLabel}{Table caption goes here.}{lcccccc}{
%    Trial & 1 & 2 & 3 \\
%    \hline
%    A & 11120739 & 11121824 & 11121825 \\
%    B & 11120739 & 11121824 & 11121825 \\
%    \hline
%}
%\ref{tableLabel}


The serial code is implemented in main.f90. The code traces $NPS$ histories and tallies the fraction of particles reflected, absorbed, and transmitted by a 1-dimensional shield of length $D=1$ cm. The shield is composed of $nmat = 100$ isotopes. Material $imat$ has arbitrary atom density
\eq{
    \rho[imat] =
    0.02 \exp(-0.001 \times imat) \ \forall \ imat \ \in [1, nmat]
}
The arbitrary energy-dependent per-atom (``microscopic'') total interaction cross section for material $imat$ is
\eq{
    \sigma_t[imat](E) =
    \sin(imat \times \ln(E)) + 1.1 \ \forall \ imat \ \in [1, nmat]
}
where $E$ is the particle energy. The microscopic scattering cross section is
\eq{
    \sigma_s[imat](E) =
    \sigma_t[imat](E) \times 0.8 \times
        (\sin(4\times imat\times\ln(E)) + 1)
    \ \forall \ imat \ \in [1, nmat]
}

The ``macroscopic'' interaction cross sections (per centimeter) are then
\eq{ \label{eq-macro-t}
    \Sigma_t(E) =
    \sum_{imat=1}^{nmat}
        \sigma_t[imat](E) \rho[imat]
}
\eq{ \label{eq-macro-s}
    \Sigma_s(E) =
    \sum_{imat=1}^{nmat}
        \sigma_s[imat](E) \rho[imat]
}

The trace step for each history is a series of ``segments.'' Each segment samples the distance the particle travels (based on $\Sigma_t(E)$) and then determines whether the particle has left the system due to escape (reflection or transmission) or absorption (based on $\Sigma_t(E)$ and $\Sigma_s(E)$). If the particle has not left the system, a new flight-direction and energy are sampled. When the particle leaves the system, the reflection, transmission, or absorption tally is incremented.

The code strides the random number sequence so that the results are exactly reproducible in MPI. Each history has 150 random numbers reserved for it. The ``advance\_rng'' subroutine calls the random number generator to exhaust the previous history's stride before the next history begins. The CPU overhead of this striding is estimated less than 5\%.

The serial code uses batched logging and retracing. For each segment, $\sigma_t[imat]$, $\sigma_s[imat]$, $\Sigma_t$, $\Sigma_s$, $dtc$ (the sampled distance to collision or distance to escape), and a loss flag are recorded. See Table~\ref{log-format}. The logs are recorded in segment-major Fortran arrays, i.e., $log\_micro(imat, iseg)$, $log\_Sigma(iseg)$, $log\_dtc(iseg)$, $log\_loss(iseg)$.

\inserttable{log-format}{Example of logged information. History 1 has 3 segments and terminates with transmission; history 2 has 1 segment and terminates with absorption; history 3 has 2 segments and terminates with reflection.}{lcccccccccccccccc}{
    $iseg$ & $\sigma$ & $\Sigma$ & $dtc$ & iloss \\
    \hline
    1 & \{$\sigma[imat=1](E_1)$, ..., $\sigma[imat=nmat](E_1)$\} & $\Sigma(E_1)$ & 0.251 & 0 \\
    2 & \{$\sigma[imat=1](E_2)$, ..., $\sigma[imat=nmat](E_2)$\} & $\Sigma(E_2)$ & 0.715 & 0 \\
    3 & \{$\sigma[imat=1](E_3)$, ..., $\sigma[imat=nmat](E_3)$\} & $\Sigma(E_3)$ & 0.212 & 2 \\
    \hline
    4 & \{$\sigma[imat=1](E_4)$, ..., $\sigma[imat=nmat](E_4)$\} & $\Sigma(E_4)$ & 0.816 & 3 \\
    \hline
    5 & \{$\sigma[imat=1](E_5)$, ..., $\sigma[imat=nmat](E_5)$\} & $\Sigma(E_5)$ & 0.005 & 0 \\
    6 & \{$\sigma[imat=1](E_6)$, ..., $\sigma[imat=nmat](E_6)$\} & $\Sigma(E_6)$ & 0.081 & 1 \\
    \hline
}

Each batch consists of $nbatch = 100$ histories. The serial calculation time is generally insensitive to $nbatch$. Indeed, batching is not necessary in the serial implementation, but is important for efficient decomposition and CUDA acceleration.

Retracing each segment involves two steps per model. First, the macroscopic cross sections in the retrace model are calculated using Eqs.~\ref{eq-macro-t}--\ref{eq-macro-s}. In the retrace model, the atom densities are slightly different than in the trace model, so the macroscopic cross sections will be slightly different. Each model corresponds to a perturbation in one material, so $nmodel = nmat$. The macroscopic cross section calculation is effectively a dot product of vectors with length $nmat$.

Second, the weight multiplier for the tally is corrected by a ratio of probability-density functions. These functions depend on whether the particle scatters ($log\_loss(iseg) = 0$), is absorbed ($log\_loss(iseg) = 3$), or escapes ($log\_loss(iseg) \in [1,2]$). Each of these calculations requires arithmetic and exponential operations. If the particle is lost, the corresponding tally is incremented and the weight multiplier is reset to unity.

The developers established a typical test problem with $NPS = 10^5$ and $nmat = nmodel = 100$. Estimates of the serial runtime dedicated to each portion of the calculation are listed in Table~\ref{serial-times}.

\inserttable{serial-times}{Estimated breakdown of a 8.2-second serial calculation with $NPS = 10^5$ and $nmat = nmodel = 100$.}{lcccccccc}{
    Step & Time [s] \\
    \hline
    Trace & 0.2 \\
    Calculate retrace $\Sigma$ & 5.5 \\
    Calculate retrace wgt & 0.8 \\
    Logging, etc. & 1.4 \\
    \hline
}
